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Abstract 

We report on a high-precision calculation of the Bhabha process in Quantum Elec- 
trodynamics, of interest for precise luminosity determination of electron-positron 
colliders involved in R measurements in the region of hadronic resonances. The 
calculation is based on the matching of exact next-to-leading order corrections 
with a Parton Shower algorithm. The accuracy of the approach is demonstrated in 
comparison with existing independent calculations and through a detailed analysis 
of the main components of theoretical uncertainty, including two-loop corrections, 
hadronic vacuum polarization and light pair contributions. The calculation is imple- 
mented in an improved version of the event generator BABAYAGA with a theoretical 
accuracy of the order of 0.1%. The generator is now available for high-precision 
simulations of the Bhabha process at flavour factories. 
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next-to-leading order corrections; two-loop corrections; Parton Shower 
pacs: 12.15.Lk, 12.20.-m, 13.40.Ks, 13.66.Jn 



1 Introduction 



The determination of the R ratio in electron-positron annihilation, denned 
as R = a(e + e~ — > hadrons)/<r(e + e~ — > fi + fi~), is a classical measurement of 
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particle physics and still a quantity of deep interest in modern research about 
the fundamental constituents of matter. The measurement of the cross section 
for electron-positron annihilation into hadrons is, in fact, an important task of 
high-luminosity colliders operating in the region of hadronic resonances, such 
as $, r-charm and B factories p]. The reason is that the R value is crucial 
for precise predictions of the hadronic contribution to (g — 2)^, the anomalous 
magnetic moment of the muon j2], and to the running of the electromagnetic 
coupling from its value at low energy up to high energies [3]. In particular, 
the QED coupling constant evaluated at Z pole olqed{.M z ) is a fundamental 
ingredient in precision tests of the electroweak theory and its uncertainty 
critically limits the bound on the Higgs boson mass in the indirect search 
through fits to precision data [I]. 

The R value has been measured by many experiments from hadron produc- 
tion threshold to the mass of the Z boson. Below the bottom threshold and, 
especially, in the energy range below 5 GeV, precision measurements of R 
are motivated by reducing the uncertainty of the hadronic contribution to 
(g — 2) M and olqed(Mz)- The experimental methods presently employed to 
determine R include direct measurements and indirect measurements of the 
hadron production cross section [5]. The first approach is followed in the re- 
cent measurements by CMD-2 and SND collaborations at VEPP-2M [6], BES 
at BEPC [7] and CLEO at CESR [8]. The indirect measurement, which makes 
use of the emission of one or more hard photons from the initial state and is 
known as radiative return [5J, is presently performed by KLOE collaboration 
at DA$NE [TO] and BABAR at PEP-II [11], and is under consideration by 
BELLE at KEK-B. The radiative return is of particular interest because it 
enables, in principle, to measure R over the full energy range, i.e. from hadron 
production threshold to the nominal centre-of-mass energy. 

Independently from the different types of error affecting the precision of the 
two methods, a common source of systematic uncertainty comes from the 
knowledge of the collider luminosity. To keep under control such an uncer- 
tainty, high-precision calculations of the QED processes e + e~ — > e + e~ , fi + yT , 
77, and relative Monte Carlo generators, are required. Among the QED pro- 
cesses, the large-angle Bhabha scattering is of particular interest because of 
its large cross section and its clean experimental signature. To simulate the 
experimentally relevant distributions and calculate the cross section of the 
Bhabha process, KLOE and CLEO collaborations make use of the QED Par- 
ton Shower generator BABAYAGA, developed in Refs. [T2"iri3] with a precision 
target of 0.5%. The Monte Carlo MCGPJ [Hj, which includes exact 0{a) correc- 
tions supplemented with leading logarithmic higher-order contributions and 
has an estimated accuracy of about 0.2%, is presently used at VEPP-2M to 
monitor the collider luminosity. To keep under control the theoretical preci- 
sion, other codes, such as the 0{a) generator of Ref. [15] (based on Ref. [16]) 
and the Yennie-Frautschi-Suura tool BHWIDE [17], are also employed by the 
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experimental collaborations. A generally good agreement between theory and 
data, as well as between the results of the different generators, is observed, 
thus confirming the precision claims of the respective calculations [TU|14lll8j . 



Nevertheless, further progress in the calculation of radiative corrections to 
QED processes and, in particular, the development of precise large-angle 
Bhabha generators are still required. This is motivated by a number of reasons. 
First, the total luminosity error quoted by KLOE is presently 0.6% [ID] , where 
the dominant source of uncertainty comes from theory, i.e. from the 0.5% phys- 
ical precision inherent the BABAYAGA generator. The reduction of such an error 
demands progress on the Bhabha theory side. Secondly, the measurement of 
the hadronic cross section in the 7r + 7r~ channel at VEPP-2M has achieved a 
total systematic error of 0.6—1% in the region of the p resonance [6], which 
requires, in turn, an assessment of the collider luminosity at the level of 0.1%. 
Last but not least, precision measurements of R trough radiative return at 
high-luminosity e + e~ storage rings KEK-B and PEP-II are already performed 
or foreseen in the near future, as previously mentioned. These facts are also 
among the motivations of the recent efforts in the direction of complete two- 
loop calculations to Bhabha scattering [19.20.21 ,22 23.24] . The need for keep- 
ing under control accurately radiative corrections to QED processes has been 
recently reinforced by the update of the e + e~ — > 7r + 7r~ cross section by SND 
collaboration at VEPP-2M. This reanalysis [25], which leads to a decrease of 
the measured cross section by two systematic errors in average with respect 
to the previous one, became necessary due to a flaw in the Monte Carlo gen- 
erators previously used in data analysis to compute radiative corrections to 
e + e~ — > 7r + 7r~ and e + e~ — > p + fi~. 

The aim of the present paper is to report on a high-precision calculation of 
radiative corrections to the Bhabha process, in order to improve the theoret- 
ical formulation of the original BABAYAGA generator down to O(0.1%). The 
approach is based on the matching of exact next-to-leading-order corrections 
with resummation through all orders of a of the leading contributions from 
multiple soft and collinear radiation, which are taken into account according 
to a QED Parton Shower. Emphasis is also put on the impact of higher- 
order non- leading- log corrections, such as a 2 L (with L collinear logarithm), 
and non-photonic light pair contributions, to arrive at a sound estimate of the 
overall theoretical error. A critical comparison of the formulation presented in 



It is worth stressing that the matching procedure here developed allows for a fully 
exclusive generation of multiple-photon radiation. With respect to different exclu- 
sive exponentiation approaches, such as the Yennie-Frautschi-Suura method [26] . 
the present approach differs in some implementation details and, in particular, in 
the resummation of non-infrared single collinear logarithms. However, since both 
formalisms coincide at first order, the differences start at 0{c?) and are not in- 
frared sensitive. 
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this paper with existing two-loop calculations allows to put on a more quan- 
titative ground the estimated theoretical uncertainty. 



The paper is organized as follows. In Sect. [2] we describe in detail the matching 
of next-to-leading-order corrections with Parton Shower, which the old version 
of BABAYAGA was based on. In Sect. [3] the predictions for large-angle Bhabha 
scattering of the improved version of BABAYAGA are compared with those of 
independent generators, both for integrated cross sections and differential dis- 
tributions of experimental interest. In Sect. H] different sources of theoretical 
uncertainties are investigated: vacuum polarization uncertainty, approximate 
treatment of two real photon emission, light pair corrections, missing virtual 
plus soft corrections to one real photon emission. Section 14.41 is devoted to the 
comparison of the formulation implemented in BABAYAGA, expanded at 0(a 2 ), 
with two existing calculations of next-to-next-to-leading-order corrections to 
Bhabha scattering: the purely photonic contribution (virtual plus soft-real 
corrections) and the two-loop Nj = 1 complete calculation in the soft-pair 
approximation for real pair production. These comparisons corroborate the 
claimed physical precision of BABAYAGA of 0(0.1%). It is worth stressing that, 
even with a complete two-loop calculation at hand, the effect of higher order 
corrections is still relevant on the scale of 0.1% accuracy, as proved in Sect. 14. ll 
Conclusions and possible developments are given in Sect. EJ 



2 Matching next-to-leading corrections with Parton Shower 

In this Section we discuss the consistent inclusion of an exact fixed order 
calculation in a cross section resumming the leading corrections up to all 
orders of perturbation theory, in a QED Parton Shower (PS) approach. The 
algorithm described below is now implemented in the new version of the event 
generator BABAYAGA [27], at present only for the Bhabha process. 

The matching of the two calculations is a non trivial task and it has to avoid 
the double counting at first order in a of the leading corrections already ac- 
counted for by the PS. 

A general expression for the cross section with the emission of an arbitrary 
number of photons, in leading-log (LL) approximation, can be cast in the 
following form: 



where U(Q 2 ,e) is the Sudakov form-factor accounting for the soft-photon (up 
to an energy equal to e in units of the incoming fermion energy E) and virtual 
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emissions, e is an infrared separator dividing soft and hard radiation and Q 2 is 
related to the energy scale of the process. \M n ,LL\ 2 is the squared amplitude in 
LL approximation describing the process with the emission of n hard photons, 
with energy larger than e in units of E. d$ n is the exact phase-space element 
of the process (divided by the incoming flux factor), with the emission of 
n additional photons with respect to the Born-like final-state configuration: 
considering Bhabha scattering and defining pi, p%, ps, Pa and ki (i = 0, • • • , n) 
as the initial-state e~ and e + , the final-state e~ and e + and photons' momenta 
respectively, d<& n = x flux reads 



1 cPp3 d 3 p 



£ 5 A (pi +P2-P3 - Pa) 



(2tt) 



3n + 2 



(ftps dPp 



if n = 

2k" 

% 

if n > 



(2) 



The cross section da^ L of Eq. (PJ is independent of the infrared separator e, 
provided it is sufficiently small. 

According to the factorization theorems of soft and/or collinear singularities, 
the squared amplitudes in LL approximation can be written in a factorized 
form. In the following, for the sake of clarity and without loss of generality, we 
write photon emission formulas as if only one external fermion radiates. We 
are aware that it is a completely unphysical case, but it allows to write more 
compact formulas, being the generalization to the real case straightforward 
when including the suited combinatorial factors. With this in mind, the one- 
photon emission squared amplitude in LL approximation can be written as 

l2 OL 1 + Z 1 2 87T 2 . . 

\M ltLL = 7T1 1 (*) ^0 ~r?2 7l V 3 



where 1 — z is the fraction of the fermion energy E carried by the photon, 
k is the photon four-momentum, I(k) is a function describing the angular 
spectrum of the photon and P(z) = (1 + z 2 )/(l — z) is the Altarelli-Parisi 
e — > e + 7 splitting function. In Eq. ([3]) we observe the factorization of the 
Born squared amplitude and that the emission factor ^P(z)I(k) E 2^n- Z ) can 
be iterated for each photon emission, up to all orders, to obtain | M. n ,LL,\ 2 - It is 
worth noticing that d 3 k/k° = (1 — z)E 2 dQ 1 dz and that in the collinear limit 
the cross section of Eq. (JT]) reduces to the cross section calculated by means 
of the QED PS algorithm described in Refs. plfT3] . 

The Sudakov form factor U(Q 2 ,e) reads explicitly 

n(QV)=exp(-^-/ + L'V L' = log^, /+= J dzP(z) (4) 
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The function I(k) has the property that / d0 7 J(A;) = \ogQ 2 /m 2 and allows 
the cancellation of the infrared logarithms. 

The cross section calculated in Eq. (pQ) has the advantage that the photonic 
corrections, in LL approximation, are resummed up to all orders of perturba- 
tion theory. On the other side, the weak point of the formula (CD) is that its 
expansion at 0(a) does not coincide with an exact 0(a) (NLO) result, being 
its LL approximation. In fact 



|^o| 2 ^o + |-Mi,ll| 2 *£i 



- ^- h log — 

= [1 + C a>LL ] \M \ 2 d<!> + \M 1 , LL \ 2 d^ 1 (5) 



whereas an exact NLO cross section can be always cast in the form 

da a = [1 + CJ \M \ 2 d<& + \M x \ 2 d^ (6) 

The coefficient C a contains the complete virtual 0(a) and the 0(a) soft- 
bremsstrahlung squared matrix elements, in units of the Born squared ampli- 
tude, and |A^i| 2 is the exact squared matrix element with the emission of one 
hard photon. We remark that C 0) ll has the same logarithmic structure as C a 
and that I-M^lzJ 2 has the same singular behaviour of |A^i| 2 . 

In order to match the LL and NLO calculations, we introduce the correction 
factors, which are by construction infrared safe and free of collinear logarithms, 



F SV = l + (C a -C a>LL ), F H = 1 ~l~ |9 (7) 





2 - \M liLL 
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\Mx, LL 


2 



and we notice that the exact 0{a) cross section can be expressed, up to terms 
of 0(a 2 ), in terms of its LL approximation as 

da a = F sv (l + C a , LL )\M \ 2 d<l> + FhIM^llI 2 **! (8) 



Driven by Eq. (jSj), Eq. ([T]) can be improved by writing the resummed cross 
section as 

oo i / n \ 

= F sv n(Q 2 , e) E - II F h, \M n , LL \ 2 d$ n (9) 

n=0 U - \i=0 I 



The correction factors Fh,% follow from the definition Eq. ([7]) for each photon 
emission. The expansion at 0(a) of Eq. (jHJ) coincides now with the exact NLO 
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cross section Eq. ([6]) and all higher order LL contributions are the same as in 
Eq. ©. 



Eq. (Q is our master formula for the matching between the exact 0(a) calcu- 
lation and the QED resummed PS cross section, according to which we also 
generate events. The extension of the matching formula Eq. (jH]) to the realistic 
case, where every charged particle radiates photons, is almost straightforward. 
We would like to remark also that the LL cross section of Eq. ([I]) is by con- 
struction positively defined in every point of the phase space, whereas the 
correction factors of Eq. ([7]) can in principle make the differential cross sec- 
tion of Eq. negative in some point, namely where the PS approximation is 
less accurate (e.g. for hard photons at large angles). Nevertheless, we verified 
that this never happens when considering typical event selection criteria for 
luminosity at flavour factories. 

It is useful to present, in the realistic case, the expression of the function I(k), 
which describes the leading behaviour of the angular spectrum of the emitted 
photons, accounting also for interference of radiation coming from different 
charged particles: 

m = £ « E * (10) 



where pi is the momentum of the external fermion I, rji is a charge factor 
equal to +1 for incoming particles or outgoing antiparticles and equal to -1 
for incoming antiparticles or outgoing particles, k is the photon momentum, 
E-y is its energy and the sum runs over all the external fermions. The function 
I(k) does not depend on the photon energy. Given Eq. ffTUl) . one can quite 
easily convince that the more convenient choice for the Sudakov form factor 
scale Q 2 is L' = log % = log — 1 = L — 1 where s, t and u are the 
Mandelstam variables of the process and m is the electron mass. 

The exact squared amplitude for the emission of a real photon |.Mi| 2 has been 
calculated by hand with the help of FORM [28J and successfully cross checked 
with Ref. [16] and with the output of the ALPHA algorithm |29j . 

The exact 0(a) soft plus virtual corrections to the Bhabha scattering have 
been taken from Ref. [30]. The soft plus virtual cross section reads 



da sv = da sv + da sv + da sv 

do a sy = dal[2(/3 + (3 mt ) logs + C>] (11) 

where i is an index for s, t and s-t subprocesses contributing to the Bhabha 
cross section, (3 = ^[log(s/m 2 ) — 1], /3 int = —\og(t/u) and the explicit ex- 
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pression for C % F can be found in Ref. [30]. We notice that in Eq. (TTTj) the terms 
coming from s, t and s-t interference diagrams are explicitly given. 

In the following, we discuss the implementation of the vacuum polarization 
effects in BABAYAGA. Some technical aspects about Eq. (its independence 
from e, the mapping of the momenta needed for n > 2 and the importance 
sampling of the final-state collinear singularities) are discussed in Appendix lAl 



2. 1 Vacuum polarization and Z exchange contributions 

Besides the photonic radiative corrections considered in the previous Sec- 
tion, also the vacuum polarization effects must be included in the master 
formula flU]), in order to reach the required theoretical accuracy. They are ac- 
counted for by replacing the fine structure constant a = a(0) with a(q 2 ) = 
Qf/(1 — Aa(q 2 )), where Aa(q 2 ) is the fermionic contribution to the photon self- 
energy: the leptonic and top-quark one-loop contributions can be calculated 
analytically in perturbation theory, while the remaining five quarks (hadronic) 
contribution, Aa^ dr , has to be extracted from data. To evaluate Aa^} dr we 
use the HADR5N routine by F. Jegerlehner [3T3~T| . 

Setting r s = a(s)/a and r t = a(t)/a (s and t are the Mandelstam invari- 
ants), we include vacuum polarization in the Born matrix element, which is 
proportional to a, by rescaling the s and t channel amplitudes, namely 

\M \ 2 = \M ,s + M ,t\ 2 -> \Mo,vp\ 2 = \M , s r s + M , t r t \ 2 (12) 

Going beyond the Born-like approximation, the cross section corrected at 0{a) 
including also vacuum polarization can be written as o~y P = aoyp + o~s V + a%, 
where a^y and a% are the soft plus virtual and the hard photon 0(a) cor- 
rections of photonic origin. We can go a step further and include vacuum 
polarization in those terms, in order to include also part of the 0(a 2 ) fac- 
torizable corrections. The hard emission matrix element is the sum of eight 
amplitudes where the real photon is attached to a s or t channel-like diagram. 
As in Eq. (fl"2l . we rescale those amplitudes by r s and r t , respectively. In or- 
der to guarantee the cancellation of the infrared separator e between ag V and 
o~x, also in ag V of Eq. (TTTj) . we rescale the s, t and s-t interference contribu- 
tions with the appropriate vacuum polarization factor, namely r 2 , r 2 and r s r t . 
Finally, the vacuum polarization improved amplitudes and cross sections are 
used as building blocks of the master formula Q). 

Furthermore, we add to the Born amplitude also the Z exchange diagrams: 
their effect is really tiny and negligible at low energies, but can become more 
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important (up to 0.1%) around 10 GeV when considering wide angular accep- 
tance regions. 



3 Numerical results 



In order to test the internal consistency of the formulation described above 
and to quantify the physical precision of the improved version of the BABAYAGA 
generator, we performed a number of Monte Carlo simulations of those exper- 
imental observables which are relevant for luminosity measurements at e + e~ 
flavour factories. To model the acceptance cuts used by the experimental col- 
laborations, we considered four different set up defined by the following selec- 
tion criteria 

(a) v/i = 1.02 GeV, E min = 0.408 GeV, 20° < 0± < 160°, ^ max = 10° 

(b) v/i = 1.02 GeV, E min = 0.408 GeV, 55° < ± < 125°, £ max = 10° 

(c) v/i = 10 GeV, E min = 4 GeV, 20° < 0± < 160°, £ max = 10° 

(d) v/i = 10 GeV, E min = 4 GeV, 55° < 0± < 125°, £ max = 10° 

where E min = 0.8 x E beam is the energy threshold for the final-state elec- 
tron/positron, 9± are the electron/positron scattering angles and £, max is the 
maximum allowed acollinearity. The set up (a) and (b) are of interest for $ 
factories, while set up (c) and (d) refer to B-factories. In both wider 
and a tighter angular acceptance are considered, in order to study the depen- 
dence of the radiative corrections from the acceptance criteria. The energy and 
acollinearity cuts are very similar to those considered in previous simulations 
and tend to single out quasi-elastic Bhabha events. It is understood that the 
energy of final-state electron/positron corresponds to a so-called "bare" event 
selection (i.e. without photon recombination), which resembles realistic data 
taking at flavour factories. 



3.1 Integrated cross sections and technical tests 



A first meaningful test of the correct matching of NLO corrections with PS is 
to prove independence of the predictions for the QED corrected cross section 
from variation of the soft-hard separator e required by the PS algorithm. 
This is successfully demonstrated in Fig. [TJ which shows the Bhabha cross 
section, obtained according to Eq. Q and in the conditions of set up (b), as a 
function of e from 10 -3 to 10~ 7 . A priori, one should expect compatibility of 
the calculated cross section against e variation at an accuracy level of 0(ae). 
This is clearly seen to be valid in Fig. [TJ when looking at the relative difference 
between the cross section predictions as e varies. Also for the 0(a) corrected 
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cross section, which is an important ingredient of the present formulation 
and a component of the following discussion, independence from e has been 
successfully checked in the limit of sufficiently small values, i.e. for e variations 
in the range 10~ 3 — 10~ 7 . 



I 







le-08 le-07 le-06 le-05 le-04 0.001 0.01 

e 

Fig. 1. QED corrected Bhabha cross section as a function of the infrared regulator 
e, according to set up (b). The error bars correspond to la Monte Carlo error. 

To quantify the overall impact of QED radiation and, in particular, to evalu- 
ate the size of QED contributions at different perturbative/precision levels, we 
show in Tab. [T]the lowest-order Bhabha cross section (without and with vac- 
uum polarization), the exact 0(a) cross section as in Eq. ([6]) and the 0(a) PS 
cross section as in Eq. (JSJ), as well as the all-order PS cross section of Eq. ([1]) 
and the matched PS cross section of Eq. (Q. The four columns correspond 
to the experimental conditions previously defined. In the cross sections of 
Tab. [U we switch off the vacuum polarization effect except in the second row, 
to better study the different sources of corrections. Interestingly, from these 
cross section values it is possible to disentangle the relative effect of various 
QED contributions, namely the contribution of vacuum polarization, of exact 
0(a) corrections, of higher-order (i.e. beyond 0(a)) leading corrections in the 
pure PS approach and in the improved PS matched with 0(a) corrections, of 
non-logarithmic terms entering the 0(a) cross section and present also in the 
improved PS algorithm and, finally, of part of the sub-leading a n L n ~ x effects. 
The above per cent corrections are shown in Tab. [2] and they can be derived 
from the cross section values of Tab. CD according to the following formulae 
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(c) 


(d) 


NLO 
u a 

a PS 

u a 

rr PS 

matched 


6855.743(1) 
6976.49(4) 
6060.07(6) 
6083.59 (6) 
6086.74(7) 
6107.57(6) 


529.4631(2) 
542.657(6) 
451.523(6) 
454.503 (6) 
455.858(5) 
458.437(4) 


71.333(1) 
74.7632(6) 
59.900(1) 
60.144(1) 
60.419(1) 
60.628(1) 


5.5026(2) 
5.85526(3) 
4.4256(2) 
4.4565(1) 
4.5046 (3) 
4.5301 (2) 



Table 1 

Bhabha cross section (in nb) according to different precision levels and for the four 
set up specified in the text. 



set up 


(a) 


(b) 


(c) 


(d) 


Svp 


1.76 


2.49 


4.81 


6.41 


S a 


-11.61 


-14.72 


-16.03 


-19.57 


&HO 


0.39 


0.82 


0.73 


1.44 


* P H S 


0.35 


0.74 


0.68 


1.34 


<w 


0.04 


0.08 


0.05 


0.10 


t-non-log 
O a 


-0.34 


-0.56 


-0.34 


-0.56 


rnon-loq 


-0.30 


-0.49 


-0.29 


-0.46 



Table 2 

Relative corrections (in per cent) to the Bhabha cross section for the four set up 
specified in the text. 

From Tab. [2] it can be seen that the vacuum polarization gives a positive cor- 
rection to the lowest-order cross section of the order of 2% at the $ factories 
and of about 5-6% at the B factories, being its contribution more important 
for a tighter angular acceptance than for a wider one. This dependence of 
vacuum polarization from the detector acceptance has to be ascribed to the 
role played by the different sub-processes contributing to the Bhabha cross 
section as the angular acceptance varies. Actually, the energy scale entering 
the logarithmic dependence of the vacuum polarization is process dependent 
and is equal to the cm. energy s for the time-like s-channel sub-process and to 
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(the absolute value of) the momentum transfer \t\ for the space- like t-channel 
contribution, being on the average \t\ <C s. While for a wide acceptance the 
t-channel contribution is largely dominating, the s-channel contribution be- 
comes more and more important as the angular acceptance decreases, thus 
explaining the trend observed for the vacuum polarization correction. The ex- 
act 0(a) corrections lower the cross section of about 15% ($ factories) and of 
about 20-25% (5-factories). Higher-order contributions of the type 0(a n L n ), 
with n > 2, introduce a positive correction around 0.5-1% at the $ factories 
and at the 1-2% level at the B factories. Therefore, multiple photon correc- 
tions are unavoidable in view of the required theoretical precision, as already 
noticed in Ref. [12]. On the other hand, also non-log O(a) corrections are nec- 
essary at a precision level of 0.1%, since their contribution is of the order of 
0.5%, almost independently from the cm. energy and with a mild dependence 
from the angular cuts. This confirms a posteriori the need for matching the 
original PS formulation of BABAYAGA with NLO corrections. The effect due to 
0(a 2 L) corrections varies from 0.05% (wide acceptance) to 0.1% (tight ac- 
ceptance). Although these contributions are only approximately kept under 
control, it can be argued, and will be shown in the following Section, that 
the infrared part of 0(a 2 L) terms is correctly reproduced by the present ap- 
proach [32] and, therefore, the size of such effects can be viewed as an estimate 
of the overall physical precision, which is conservatively close to 0.1%. From 
Tab. [2] it can be also seen that the matching of NLO corrections with PS does 
not alter at the level of 0.1% the size of higher-order and NLO contributions, 
thus preserving correctly the impact of these partial effects. This conclusion 
can be inferred by comparing 5 ho with S^q and S™ on ~ l ° 9 with ^ n "' 09 , respec- 
tively. A common feature observed for all the relative corrections, with the 
exception of non-log effects, is the fact that, at a fixed cm. energy, they are 
larger (by about a factor of two) for a tighter acceptance with respect to a 
wider one. This can be understood as follows. As discussed above, a natural 
choice for the Q 2 dependence of the collinear logarithm L = log Q 2 /m 2 is 
Q 2 = st/u. This implies that at large scattering angles the collinear logarithm 
tends to log s/m 2 , while for a wide angular acceptance it can be on the average 
approximated by log|i|/m 2 , thus explaining the observed angular-dependent 
behaviour. 



3.2 Differential distributions 

Besides the integrated cross section, various differential distributions are used 
by the experimental collaborations to monitor the collider luminosity. We 
show, in Fig. [2] and Fig. [3j two distributions which are particularly sensitive 
to the details of photon radiation, i.e. the e + e~ acollinearity and the invariant 
mass distribution, in order to quantify the differences between the previous 
version of BABAYAGA (denoted as OLD in the figures) and the improved one 
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presented here (denoted as NEW). As a reference, the distributions obtained 
according to the exact 0{a) calculation are also shown. The results refer to 
set up (b). 
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Fig. 2. Acollinearity distribution according to the PS matched with 0(a) corrections 
(Eq. @, solid line), the LL PS algorithm (Eq. ([!]), dashed line) and the exact 
0(a) calculation (dotted line). The inset shows the relative differences between the 
predictions of the improved and the LL PS. Selection criteria of set up (b) are 
considered. 
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Fig. 3. The same as Fig. [2] for the e + e invariant mass distribution. 

From Fig. [2] and Fig. [3] it can be clearly seen that multiple photon corrections 
introduce significant deviations with respect to an 0(a) simulation, especially 
in the hard tails of the distributions, where they amount to several per cent. 
To make more visible the improvements introduced by the matching proce- 
dure discussed in Sect. I2J the inset shows the relative differences between the 
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predictions of the improved and original version of BABAYAGA. These differ- 
ences mainly come from non-log 0(a) contributions and, to a smaller extent, 
from 0(a 2 L) terms. Their effect is flat and at level of 0.5% for the acollinear- 
ity distribution, while they reach the some per cent level in the hard tail of 
the invariant mass distribution. As a whole, these results demonstrate that 
exact 0(a) and higher-order corrections need to be simultaneously taken into 
account for precision luminosity studies. 

3.3 Tuned comparisons 

An important step towards the estimate of the theoretical accuracy of the for- 
mulation is the tuned comparison of the improved version of BABAYAGA with 
independent precision calculations of the Bhabha process. To this end, we first 
compared the results for the integrated cross section as obtained by the im- 
proved version of BABAYAGA with the corresponding predictions of LABSPV [33] 
and BHWIDE [17J, which both rely on different theoretical ingredients. For the 
sake of comparison, the contribution of vacuum polarization has been switched 
off, in order to test just the implementation of pure QED corrections. The com- 
parison is shown in Tab. [3j for both set up (a) and (b). As can be seen, the 
predictions of the three calculations agree within 0.1%. Although we didn't 
performed detailed tests in comparison with the recently developed generator 
MCGPJ [Hj, we expect a level of agreement with that calculation at a similar 
precision level, on the basis of the comparisons between BHWIDE and MCGPJ 
reported in Ref. [H]. 



set up 


(a) 


(b) 


^BHWIDE 

Olabspv 

^matched 
"BABAYAGA 


6086.3(2) 
6088.5 (3) 
6086.61 (2) 


455.73(1) 
456.19(1) 
455.853 (4) 



Table 3 

Comparison between the predictions for the Bhabha cross section (in nb) as obtained 
with BHWIDE, LABSPV and the present version of BABAYAGA. 

Comparisons between BABAYAGA and BHWIDE at the level of differential distri- 
butions are given in Fig. H] and Fig. [51 where the inset shows the relative devi- 
ations between the predictions of the two codes, with reference to set up (b) . 
As can be seen, there is a very good agreement between the two generators, as 
the predicted distributions appear, at a first sight, almost indistinguishable. 
Looking in more detail, there is a relative difference of a few per mille for the 
acollinearity distribution (Fig. H]) and of a few per cent for the invariant mass 
(Fig. [5]), but only in the hard tails, which little contribute to the integrated 
cross section. In fact, these differences on differential distributions translate 
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into an agreement on the cross section values well below the 0.1% level, as 
shown in Tab. [3j 

On the ground of these results, we can conclude that the calculation of QED 
corrections in the two generators, although based on different approaches, 
numerically agrees very well, at the level of the required precision for accurate 
luminosity measurements at flavour factories. 
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Fig. 4. Comparison between the BHWIDE generator (dashed line) and the present 
version of BABAYAGA (solid line) for the acollinearity distribution. The inset shows the 
relative difference between the predictions of the two generators. Selection criteria 
of set up (b) are considered. 
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Fig. 5. The same as Fig. [J] for the e + e invariant mass distribution. 
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4 Estimate of the theoretical accuracy 



Since different implementations of radiative corrections beyond exact 0(a) 
contributions differ by higher order effects, the results of the comparisons 
quoted in the previous section between BABAYAGA and other event generators 
give a hint of the missing radiative corrections which dominate the theoretical 
accuracy. From the investigations of the previous Section, related to the higher 
order terms inherent in the theoretical formulation of BABAYAGA and to tuned 
comparisons with independent codes, an estimate of the physical precision in 
the calculation of the radiative corrections of the order of 0.1% can be inferred. 

This guess can be verified by a comparison with a complete two-loop calcula- 
tion, since all neglected terms in the formulation described in Sect. [2] start at 
two-loop order. Being a full two-loop calculation not presently available, the 
aim of the present Section is to try to estimate the impact of the uncertainties 
at order a 2 within the realistic set up considered in this paper, by comparing 
with some of the available calculations in the literature. 

Another important source of error is the uncertainty on the hadronic contribu- 
tion to the vacuum polarization Aa^ dr , which will be quantified in Sect. 14.71 
The origin of this error is intrinsically non-perturbative because Aa^ dr can 
not be calculated perturbatively around the hadronic resonances and must 
be calculated via dispersion relations by means of data [3|. As discussed in 
Sec. I4.7[ this error will be the dominant one close to the J/^> resonances. 

We start by considering the theoretical error of perturbative origin, at 0(a 2 ). 
The pure 0(a 2 ) content of our master formula (jHJ) can be read by expanding 
it and can be cast in the following form 

9 9 9 9 

a =cr sv + a SVH + a HH (16) 



where a% v contains all the a 2 virtual and real contributions without photons 

2 

with energy larger than e, <Jsvh contains all the virtual and real contributions 
with at least one photon with energy larger than e and ojj H is the contribution 
with two real photons with energy larger than e. In BABAYAGA, each of the 
three terms is affected in principle by an error and in the following we try to 
study the impact of this error on the integrated cross section to establish the 
theoretical accuracy of the master formula (151) . 
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4-1 Corrections beyond two-loop 



Before starting to discuss the theoretical error, we would like to prove that 
the LL radiative corrections beyond a 2 are still very important, at least when 
considering differential distributions. This means that even if the complete 
two-loop perturbative calculation will be fully available, it would be desirable 
a matching with the resummation of all the remaining LL corrections. 

20 




10 



£ (cleg.) 



Fig. 6. Impact of a 2 (solid line) and higher order (dotted) corrections on the 
acollinearity distribution, in set up (b). 

This can be demonstrated by comparing the fully resummed cross section of 
Eq. (jnj) with its expansion up to order a 2 , which, relying on the results of the 
next Sections, can be considered a really good approximation of the complete 
two- loop calculation. As in Sect. 13.31 here we switch off vacuum polarization 
effects. 



In Fig. [6j the effect of the higher order corrections, dominated by a 3 con- 
tributions, is shown in comparison with that of the a 2 corrections on the 
acollinearity distribution for set up (b): as can be seen, the a 3 effect can be as 
large as 10% in the phase space region of soft photons emission, corresponding 
to small acollinearity angles with almost back-to-back final-state fermions. 



4-2 Error on two real photons emission cross section 



The cross section cr^- ^ for the emission of two real photons in BABAYAGA is 
approximated by the integral over the phase space of F^hF2,h \M.2,Lh\ 2 d<&2- 
By means of the ALPHA [29J algorithm, the two real photon amplitude can be 
calculated exactly and then integrated over the (exact) phase space c?$ 2 . In 
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Tab. HI the two real photons cross section as obtained with the approximation 
of Eq. (J9j) and with the exact matrix element are compared. We choose e = 
1 x 10 -4 and we successfully verified that the difference of the cross sections 
A is independent of the e choice, as expected. The error 8^ r H is measured in 
units of the Born cross section. 



set up 


<h (nb) 


ct alpha ( nb ) 


A (nb) 


8% H (%) 


(a) 


2189.37(2) 


2189.4(2) 


-0.3(2) 


-4(3) xlO" 3 


(b) 


229.197(2) 


229.20(2) 


0.00(2) 


0(4) xl0~ 4 


(c) 


44.0719(5) 


44.072(5) 


-0.00(5) 


0(7) xl0~ 3 


(d) 


4.22839(4) 


4.2286(5) 


0.0002(5) 


4(9) xl0~ 3 



Table 4 

Error on two real photon emission cross section. 

In Tab. HI we remark the small (and negligible) error, showing a posteriori 
that having corrected the real LL emissions with the factors Fn,i in the master 
formula gives an extremely good approximation of the exact 0(a 2 ) squared 
amplitude, at least within the considered event selection criteria. 

4-3 Light pair corrections 

An 0(a 2 ) contribution which is not included in Eq. ffl3l) . even approximately, 
are the so-called pair corrections. Starting from the photonic one loop cor- 
rected diagrams, the virtual pair diagrams can be obtained by inserting a 
fermion loop in the correcting photon propagator. In order to estimate the 
size of the virtual pair corrections (VPC), we use the formulae of Ref. [31] 
(approximated at LL accuracy), for t channel Bhabha. VPC develop terms 
of the order of a 2 L 3 which are then cancelled when also the real pair emis- 
sions (RPC) are included. We also estimate the size of the real pairs using 
the soft approximation of Ref. [351136] and setting the maximum pair energy 
AE = 0.2 x y/s, which is compatible with the requirements of set up (a)- 
(d) for the energy of the final-state leptons. In Tab. El the virtual and real 
pair cross section are reported, and the errors 5^[ irs due to missing them in 
BABAYAGA is measured in units of the Born cross section. It is worth noting 
that the relative contribution of light pairs emission is at the level of a few 
0.01%. 

Here we consider only electron pairs. The contribution of muon pairs and 
hadronic pairs is expected to be one order of magnitude smaller with respect 
to electron pair production. This problem was investigated in Ref. [37] for 
the small angle Bhabha scattering at LEP1, where the typical momentum 
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set up 


2 

cr?>n^ fnb) 
VPC \ J 


2 

u RFC \ J ) 


5 err (%) 

w pairs \ /v J 


(a) 


-4.bUo(o J 


o.oUo^o J 


-u.uiy 


(b) 


-0.5698(1) 


0.4375(1) 


-0.025 


(c) 


-0.1385(1) 


0.1154(1) 


-0.032 


(d) 


-0.01542(1) 


0.01320(1) 


-0.040 



Table 5 

Error due to the virtual and real electron pair corrections. 

transfer is of the order of the energies involved at flavour factories. In Ref. [37] 
the global effect of muon and hadronic pairs was safely estimated at the level 
of 30% of the contribution of electronic pairs. 

4-4 Comparisons with virtual plus soft two-loop calculations 

In order to establish the error on crgy of Eq. (fl~3|) . we compare it with recent 
calculations appeared in the literature. As mentioned in the introduction, there 
has been important progress towards the calculation of the full QED two-loop 
corrections (NNLO) to Bhabha scattering. An exhaustive report of the sta- 
tus of the two-loop QED corrections to Bhabha scattering can be found in 
Ref. [19] . What is actually available is the complete two- loop virtual photonic 
correction in the approximation of neglecting 0(m 2 /Q 2 ), where Q 2 stands 
for one of the Mandelstam invariants s, t and u [21]. The real radiation con- 
tribution is treated in the soft photon approximation and the calculation is 
differential in the electron scattering angle. The results have been confirmed 
by independent calculations [20,23j. In Ref. [23] the electron mass terms have 
been included but the two-loop box diagrams have been neglected. Work is 
in progress towards the calculation of massive two-loop box diagrams |24j . 
Another ingredient towards the two-loop complete calculation is the virtual 
fermionic contribution Nf = 1 [22] including all finite fermion mass effects. In 
Ref. [23J also real pair production [35J and photon emission in the soft limit 
have been introduced verifying the cancellation of infrared singularities as well 
as terms proportional to L 3 between virtual and real corrections. 

We remark that, in order to achieve a NNLO accuracy in Bhabha scattering, 
also the exact two real photons corrections and the exact one-loop corrections 
to the one real photon emission process must be included [36 38J. The complete 
calculation of the latter in particular is still missing in the literature for the 
large-angle Bhabha scattering. 

The cross section o~gy of Ecj. (1131) can be directly compared with the results 
of Ref. [21~f 221123) . in order to quantify the size of the missing terms. In the 
following subsections the uncertainties inherent to the classes of pure photonic 
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and Nj = 1 corrections are separately investigated numerically. 

As a first step we write explicitly cr^y, which is derived from the first term 
(to = 0) of the infinite sum in Eq. (Q. In order to show the s, t and interference 
contributions, we define 



da sv = dcr s,sv + da t,sv + da st,sv = ( E * + E t + E st )da 
dcr = dcr sfi + dcr tfl + do stfl = (B s + B t + B st )dcr (14) 

Truncating every factor in Eq. (Q, improved with vacuum polarization effects 
as described in Sect. 12.11 at 0(a 2 ) we get 



da 



1 + V 




+ (E s - VB s )r 2 s + {E t - VB t )r 2 t + {E st - VB st )r s r t 
(B s r 2 s + B t r 2 t + B st r s r t 



(15) 



where V = —(2a/ir)I + L' is the 0(a) term of the Sudakov form factor, £7, and 
Bi have been defined above and r s i are the vacuum polarization corrections 
(including only electron loop to be consistently compared with the Nf = 1 
results) for the s and t channels. If we define 1/(1 — Aa(q 2 )) = 1/(1 — 5 q 2), 



s- 



rf and r s r t read: 



r 2 = 1 + 2S S + 3S 2 
r 2 = l + 25 t + 35 2 

r s r t = l + 5 s + 5 t + 5 2 + 5 2 + 5 s 5 t (16) 
We would like to stress that Si are calculated at one-loop order. 
Retaining only terms up to O (a 2 ), Eq. ( fi~5{) reads 



dcrsv _ 1 
da 

+ V+(E S - VB S ) + (E t - VB t ) + (E st - VB st ) 
+ 2{B s S s + B t S t ) + B st {S s + S t ) 
+ 1/2V 2 

+ {E a - VB S )S S + {E t - VB t )S t + {Est - VB st ){S s + S t ) 

+ 3{B S S 2 + B t S 2 ) + B st {S 2 + S 2 + S s S t ) 

+ V[{E S - VB S ) + {E t - VB t ) + {E^ - VBst] 

+ V[2{B s S s + BtSt)+Bst{Ss + St)} 

+ [{E s - VB S ) + {E t - VB t ) + {E^ ~ VB st ] x 
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x [2(B S S S + B t 5 t ) + B st (S s + S t )} 



(17) 



The first line of the previous equation is the Born contribution, the second 
line is the photonic one loop soft plus virtual correction (notice that it is equal 
to E s + E t + E st because B s + B t + B st = 1), the third line is the vacuum 
polarization correction at O(a) and the remaining lines represent o~g V . 



4-4-1 Two-loop photonic corrections 

After switching off the terms coming from vacuum polarization contributions 
(5 S — S t — 0), the pure 0(a 2 ) term of the above equation can be compared with 
the analytical spectrum of Ref. [2T]. In fact the approximation s, t, u » m 2 
is fulfilled for the event selections (a), (b), (c) and (d) considered in Sect. El 
Since all infrared terms are factorized, all differences between Eq. (fTTj) and the 
calculation of Ref. [21] are not expected to be infrared sensitive, apart from 
spurious terms in Eq. (fTTl) suppressed by coefficients of the order of m 2 /Q 2 \^\. 
The infrared behaviour of the difference between Eq. (fTTj) and the results of 
Ref. [2T] can be seen in Fig. [7| (upper curve), where the infrared regulator has 
been allowed to scan over a range of ten orders of magnitude for set up (a). 
The dashed curve fitting the circles has the following expression: 
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Fig. 7. Difference between the corrected cross section of Eq. (|17p and the exact two 
loop photonic and Nj = 1 corrections, as a function of the infrared regulator e, 
according to set up (a). 



2 In BABAYAGA, the terms proportional to m 2 /Q 2 are only partially accounted for: 
for example, the fully massive kinematics is always considered, while the non- 
infrared 0(m 2 /Q 2 ) mass terms are neglected in the virtual one loop contributions. 
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a 



or, phot. 
SV 



— a 



Ref. Mi 



a = -4.02 ±0.01 

b = -6.7 ±0.4 

c= ±1.75447 ±0.00001 



?77 ^ Tfl^ \ 

a — L 2 log 2 £ ± b — L 2 log e + cL) a 



with cr = 6855.7 nb and L = log(s/m 2 ) = 15.2. Equation (IT81) clearly demon- 
strates that in the difference infrared sensitive contributions survive, but are 
suppressed by the factor m 2 /s. It also shows that, concerning photonic cor- 
rections, the error of BABAYAGA starts at the level of the a 2 L corrections, not 
enhanced by any infrared logarithm 



In order to numerically check that the term cL in Eq. fllSp is a true single 
collinear logarithm, Fig. [S] shows a scan of the difference between the QED 
corrected Bhabha cross section of Eq. (|T7|) and Ref. [21J (upper curve) as a 
function of the electron mass, whose values are allowed to span a range of 
eight orders of magnitude. The infrared regulator e has been fixed to 1 x 10~ 5 . 
The dashed curve fitting the circles has the following expression: 
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Fig. 8. Difference between the corrected cross section of Eq. (|17j) and the exact two 
loop photonic and Nj = 1 corrections, as a function of the electron mass, with the 
infrared regulator e = 1 X 1CP 5 , according to set up (a). 



VSV ^ - ^Ref. m = (-) ( fl lQ S (^2) + ^ a ° 

a = ±2.656 ±0.001 

b = -1.391 ±0.003 (19) 
showing that the missing two-loop photonic contributions in BABAYAGA are 
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really of the order of a 2 L. The impact on the integrated cross section within 
the realistic set up of these terms will be shown in Tab. [61 



4-4-2 Two-loop Nj = 1 corrections 

The complete Nf = 1 virtual calculation of Ref. [22] refers to three types of di- 
agrams: 1) diagrams with a fermionic loop in the tree- level photon propagator 
with a further photonic line within the loop; 2) diagrams with a fermionic loop 
in the photon propagators of the vertex and box one-loop diagrams; 3) dia- 
grams with a fermionic loop in the tree-level diagrams with an additional 
photonic vertex correction. In the literature, type 2) diagrams are referred to 
as virtual pair corrections, while type 3) diagrams can be thought of as factor- 
ized in vacuum polarization times one-loop photonic correction contributions. 

2 

The part of a^ v to be compared with two-loop Nf = 1 corrections is extracted 
from Eq. (|17|) by retaining all terms of photonic and vacuum polarization ori- 
gin. In order to compare with the results of Ref. [23], we subtract from them 
the virtual pair corrections to the t-channel Bhabha scattering already consid- 
ered in Sect. 14.31 (|34j). Furthermore, the soft real pair corrections included in 
Ref. [23J have been switched off, in order to perform a consistent comparison. 

We expect that all infrared sensitive terms have the correct coefficient, owing 
to the factorization of infrared corrections, as can be seen in Fig. [7| (lower 
curve). The analytical expression of the dotted curve fitting the diamonds is 
given by 



a 



A' 



sv 



} 



a Ref. [23] a VPC 



— \ la — L log e + oL I cr 



a = -8.72 ±0.07 

b = -0.955489 ± 0.000006 



(20) 



showing that the two-loop Np = 1 infrared structure in BABAYAGA is under 
control. 



Concerning the dependence from the (logarithm of the) electron mass, apart 
from spurious coefficients of order m 2 /Q 2 and a small residual a 2 L 3 depen- 
dence due to VPC to s-channel Bhabha, the terms containing collinear loga- 
rithms have two sources, namely squared vacuum polarization and interference 
between vacuum polarization and photonic corrections. Thus the difference 
between a^y f ~ l and the calculation of Ref. [23] (subtracted of the VPC of 
Sect. I4.3P should start from terms of the order of a 2 L 2 . In fact, the dotted 
curve fitting the diamonds of Fig. [8] has the following form: 
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a 



a\N f 



sv 



a 



Ref. 1231 



°VPC 



a 



71 



- a log d — ± b log' — ± c log — ± oV <7 



a = +0.00728 ±0.00007 
b = -0.504 ±0.005 
c = ±5.9 ±0.1 
d= -1.36 ±0.08 



(21) 



In the next Section the size of the error induced by missing exact photonic 
and N- 



f 



1 corrections in a^ v is discussed. 



4-5 Error on purely virtual plus soft two-loop cross section 



The numerical effect on the cross section, for the four different set up con- 
sidered in this study, is shown in Tab. [61 with the infrared separator e set 
to 10~ 5 . In Tab. El the A are the difference of the exact and the BABAYAGA 
a 2 virtual plus soft cross section and 5 err are their size in unit of the Born 
cross section. It is worth noticing the accidental opposite sign between S p h t 



set up 


Aphot (nb) 


A Nf= i (nb) 




W t =i (%) 


6 p r hot+N f =l ( % ) 


(a) 


0.9855 


-0.537(2) 


0.014 


-0.0078 


0.0062 


(b) 


0.1218 


-0.0442(3) 


0.023 


-0.0083 


0.0147 


(c) 


0.0149 


-0.01436(4) 


0.021 


-0.020 


0.001 


(d) 


0.0016 


-0.00099(1) 


0.029 


-0.018 


0.011 



Table 6 



Cross section difference (&- p hot,N f =i) as obtained with <jg V and the exact two loop 
corrections, for purely photonic contributions and Nj = 1. 8 p hot,Nt=l are the differ- 
ences of the cross sections in units of the Born one. 

and 5/v /= i for all event selections considered, allowing to reduce the error of 
BABAYAGA predictions when summing up everything together. 



4-6 Error on the virtual plus soft corrections to one photon real emission 



The virtual plus soft cross section with one hard photon Cg VH of Eq. ([TBI 
is obtained from the 0(a 2 ) content of the term n — 1 in Eq. Its exact 
expression for large angle Bhabha scattering is not known in the literature, 
but it was calculated for the small angle Bhabha process [38] or for the s- 
channel |39j at large angles. Relying on the LEP experience and on the purely 
virtual plus soft results described in Sects. 14.4.11 and 14.4.21 the error of the 
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BABAYAGA approximation is at the level of the infrared-safe a 2 L terms, the 
size of which can be safely taken to be smaller than the 0.05% for all the set 
up. 

4-7 Non-perturbative error induced by hadronic contribution to vacuum po- 
larization 

By means of the analysis presented in the previous Sections, we can deduce 
that the missing 0(a 2 ) (perturbative) contributions in BABAYAGA do not exceed 
the 0.1%. 

Nevertheless, besides the missing 0(a 2 ) contributions considered before, an- 
other source of theoretical error, coming from the hadronic contribution to the 
vacuum polarization, has to be carefully considered. This uncertainty has an 
intrinsically non-perturbative origin. Since the routine HADR5N, by means of 
which we calculate Aa^} dr (q 2 ) , returns also an error S^adr 011 its value, we esti- 
mate the induced error by computing the cross section with Aa^ dr (q 2 )±5had r 
and taking the difference as the theoretical uncertainty due to the hadronic 
contribution to vacuum polarization. In Tab. [7J the uncertainties A are calcu- 
lated for the event selection criteria (a)-(d) on the Born and the matched cross 
section. The quantities dyp are the errors in units of the Born cross section. 
Differently from previous investigations of vacuum polarization uncertainty 
in small-angle Bhabha scattering at LEP |40j, we consider here only the er- 
ror induced by the parameterization of hadronic loops of Refs. [31I3T] , because 
of the absence of other results able to keep under control appropriately the 
contribution of hadronic resonances for time-like momenta circulating in the 
photon self-energy. 



set up 


A B orn (nb) 


A f u ii (nb) 


S VP,Bom ( % ) 


^PJull (%) 


(a) 


-0.48 


-0.50 


-0.007 


-0.007 


(b) 


-0.00070 


-0.0014 


0. 


0. 


(c) 


0.017 


0.014 


0.024 


0.020 


(d) 


0.0033 


0.0024 


0.060 


0.044 



Table 7 

Error induced by the uncertainty on Aa^ dr (q 2 ). 

It must be remarked that the error on Aa^ dr (q 2 ) for q 2 > strongly varies 
and increases passing through hadronic resonances. This yields a worsening 
of the theoretical accuracy of Bhabha cross section, in particular close to the 
cc states in the region between 3 and 4.5 GeV which have a non-negligible 
branching ratio into electrons. Relying on the output of the HADR5N routine 
(which here prints a warning about its reliability), we verified that the error 
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on the Bhabha cross section induced by hadronic loop may reach the 0.5% in 
this region and it is therefore the limiting factor for the theoretical accuracy. 



4-8 Summary of the theoretical errors 

The size of the a 2 contributions missing in BABAYAGA and the error induced 
by the hadronic contribution to the vacuum polarization, within the realistic 
event selection criteria considered in this paper, are summarized in Tab. [HI 



\5 err \ (%) 


(a) 


(b) 


(c) 


(d) 


1 X eTr 1 
\°VP\ 


0.01 


0.00 


0.02 


0.04 


Impairs 1 


0.02 


0.03 


0.03 


0.04 


\°H,H\ 


0.00 


0.00 


0.00 


0.00 


1 X err 1 


0.01 


0.01 


0.00 


0.01 


1 X.err 1 
l°SV,-H"l 


0.05 


0.05 


0.05 


0.05 


1 X err 1 
\°total\ 


0.09 


0.09 


0.10 


0.14 



Table 8 

Summary of different source of theoretical error in BABAYAGA. 



The size of the virtual and real pair corrections and of the missing virtual 
corrections to the real photon emission process is only guessed with a safe 
estimate of their impact. We also remark that the somehow large error coming 
from vacuum polarization uncertainty in set up (c) and (d) is due to the fact 
that we consider energies around the T resonance and the data-driven routine 
HADR5N produces here larger errors. 

From Tab. [HI by summing up the absolute value of all the uncertainties, we 
can deduce that a safe estimate of the total theoretical accuracy of the new 
version of BABAYAGA (based on the master formula ([9]) and the HADR5N routine 
for the hadronic contribution to the vacuum polarization) for the calculation 
of the Bhabha cross section is 0.1% at $ factories and below 0.2% at the B 
factories. 

We remark that concerning the perturbative 0(a 2 ) error, our formulation is 
able to reach the 0.1% accuracy in all the considered set up and that, as 
previously emphasized, the accuracy can be worsened to the 0.5% level in 
proximity of the J/\l/ resonances because of the error on the non-perturbative 
hadronic contribution to the vacuum polarization. 
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5 Conclusions 



Bhabha scattering is a crucial process for precise luminosity monitoring of 
e + e~ colliders. Ongoing and future experiments at e + e~ accelerators operat- 
ing in the region of hadronic resonances, such as $, r-charm and B factories, 
require a precise luminosity knowledge, especially in view of improved mea- 
surements of the hadronic cross section in the energy region going from the 
pion pair production threshold up to 12 GeV. 

With this motivation in mind, we have presented in this paper a high-precision 
calculation of the Bhabha process in QED, with a precision target of the or- 
der of 0.1%. This theoretical accuracy is achieved by the matching of exact 
NLO corrections with higher-order leading logarithmic contributions, which 
are kept under control through all orders of a by means of a QED Parton 
Shower approach. The matching of matrix element corrections with Parton 
Showers is a topic of QCD calculations and recently different successful so- 
lutions have been proposed, opening the way towards precision calculations 
of high-energy QCD processes [U]. Our matching algorithm represents, to 
the best of our knowledge, the first example of such an application in QED. 
The theoretical precision of the approach has been demonstrated through de- 
tailed comparisons with the predictions of precise independent generators and, 
noticeably, with the results of recent calculations of two-loop corrections to 
Bhabha process, which are the frontier of modern perturbative QED. Other 
components of the theoretical luminosity error, coming from hadronic vacuum 
polarization and light pair emission, has been carefully estimated as well, in 
order to arrive at a sound and robust total error budget. In particular, the 
perturbative contribution to the error coming from missing 0(a 2 ) corrections 
has been shown to be under control at the 0.1% level for all the experimental 
selection criteria here considered, realistic for data analysis at flavour facto- 
ries. Furthermore, the non-perturbative error due to the uncertainty on the 
hadronic contribution to vacuum polarization, which critically depends on the 
centre-of-mass energy, turns out to be below the 0.05% level at $ and B facto- 
ries, whereas can increase up to the 0.5% in proximity of the cc bound states. 
The reduction of such an uncertainty, if needed, requires progress both on the 
theory and the experimental side. 

The theoretical precision for luminosity monitoring through large-angle Bhabha 
scattering at flavour factories is now comparable with the accuracy reached 
at the time of the LEP workshop working groups 95/96 for the theoreti- 
cal error underlying the high-precision small-angle Bhabha measurement at 
LEP [30132] . 

The new formulation has been implemented in an improved version of the 
event generator BABAYAGA [27] , which is available for high-precision luminosity 
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measurements at flavour factories. For example, this will allow to reduce the 
theoretical error in the luminosity measurement at the $ factory DA$NE 
from the present 0.5% to 0.1%, thus roughly halving the total luminosity 
error quoted by KLOE collaboration at DA$NE and, more generally, paving 
the road to more precise measurements of the Bhabha process at other e + e~ 
collider, such as BEPC, CESR, KEK-B, PEP-II and VEPP-2M. 

As far as possible future applications of the approach here presented are con- 
cerned, it would certainly be of interest for future data analysis extending 
the present phenomenological analysis to the production processes of muon 
and photon pairs in e + e~ annihilation, as well as to their exclusive signatures 
such as e + e~ — > /i + /i~7, because also these QED reactions are employed at 
flavour factories for precise luminosity studies [1]. Other interesting perspec- 
tives would be the applying the matching procedure to the Standard Model 
processes presently under consideration for luminosity monitoring of a future 
e + e~ collider at the TeV scale and extending the matching to incorporate also 
the full two-loop corrections. 

These developments are by now under consideration. 
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A Technical details 



Here we discuss some technical details about Eq. (jUj) and its integration, 
namely the independence from e, the mapping of the momenta needed for 
n > 2 and the importance sampling of the final-state collinear singularities. 
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A.l Independence of the master formula from the infrared separator 



Considering Eq. ([!]), its independence from e can be demonstrated analytically 
if we neglect the variation of |A^o| 2 from the center of mass energy. In this 
case, integrating over all the phase space, we get 

^exp(--/ + L')E f L*q.a = * (A.l) 



Considering now Eq. (Q, we can write the integral of the radiation factor for 
each emitted photon as ^L'(I + + c), where c is a constant not singular as e 
goes to zero and coming from the correction factor Fjj. Thus, the total cross 
section can be written as 



a 00 / n \ n 

a-Fsvexpi—I+L'^f^-) L' n (/+ + c) n a 

27r n=0 ^ 2nJ 

OL 

= F sv exp(— c)<x (A.2) 

Z7T 



which does not depend on e. It is worth noticing that having corrected each 
hard photon emission in Eq. with the factors Fgj is crucial for the good 
infrared behaviour of the integrated cross section. In Fig. [T] of Sect. [3] the 
numerical independence from the infrared separator has been shown. 



A.2 Mapping of the momenta for n > 2 



In the master formula (jSJ), the sum over all possible photon multiplicities 
is required, but the building blocks of the matched cross section (i.e. the 
squared matrix elements) are strictly defined only for or 1 photon in the final 
state. It is therefore mandatory to devise an algorithm to map a n photons 
momenta configuration to a or 1 photon configuration in order to consistently 
calculate the squared amplitudes |7W 1 2 and |.Mi| 2 . The mapping algorithm 
is not unique, but the choice among different mappings gives a higher order 
effect which has been verified to be negligible. 

When a photon multiplicity n is selected, the first step to define a zero-photon 
kinematics configuration is to associate ni photons as emitted by the initial 
state and np by the final state. The association is done according to an al- 
gorithm finding whether the photon is nearer to an initial-state or final-state 
fermion, in terms of its angle with respect to the fermion. If Kj = Y^jL\ kj 
and q = p\ + p 2 — Kj, two initial-state mapped momenta are defined such as 
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(pi,M + P2,m) 2 = I 2 and Pi m = v\ m = m2 - The mapped final-state momenta 
are denned such as (p 3)M + p i>M ) 2 = (Pi,m +P2,m) 2 > P%,m = p\m = m2 and 
p 3t M is directed along the direction of p 3 boosted in the frame where p 3 + p 4 
is at rest. The mapped momenta, which satisfy by construction momentum 
conservation and on-shell relations, are used to calculate the Lorentz-invariant 
Born squared amplitude needed to compute \M. n ,LL\ 2 in Eq. Q). 

In order to calculate the correction factors Fhj, a mapping to 1 photon con- 
figuration is needed. Suppose we are calculating Fn,i for the Z-th photon: in 
order to get the mapped momenta, we do as if only this photon is present. We 
keep pi^M = Pi, P2,m = P2 and k^ M — ki, we boost p 3 and p 4 where p 3 +p 4 is at 
rest, we calculate p 3)M and p 4)M such as (p 3 ,m + P4,m) 2 = (Pi,m +P2,m - h,M) 2 , 
p\m = p\m = m2 and p 3j M is directed along p 3 in this frame and finally we 
boost p 3i M and p 4) M back in the original frame. We obtain a set of mapped 
momenta satisfying momentum conservation and on-shell relations which can 
be used to calculate |.Mi| 2 and I-M^lzJ 2 and the correction factors F Hji . 



A.3 Importance sampling of final- state collinear singularities 

The integral over the phase space dQ n is done by means of Monte Carlo tech- 
niques. For n photons in the final state, the number of independent integration 
variables is 3n + 2, and we choose them to be the azimuthal and the polar 
angles of one of the final-state fermions and the energies and the angles of the 
n emitted photons. With this choice, <i$ n = <i$ n x flux can be written as 

where K = Y^7=i ^« and the final-state momenta are easily calculated with 
the independent variables by requiring momentum conservation and on-shell 
relations. In Eq. (IA.3I) . the angles of the final-state electron have been chosen, 
but the positron could have been equivalently chosen. 

The generation of the independent variables is done according to the peaking 
structure of the function to be integrated. Here we discuss in some detail the 
sampling of the final-state collinear singularities. The peaks of the differential 
cross section can be read from | M. n ,LL\ 2 d^n (the correction factors F Hji tend 
to 1 in the singular regions), namely 

\M niLL \ 2 al<S> n <xdk%--- dk^dQ^i ■ ■ ■ dVt in dVt 3 1 1(h) ■ ■ ■ I(k n ) (A.4) 

fc l ' ' ' fc n 



30 



In the previous equation the infrared and collinear singularities for each emit- 
ted photon are evident. If infrared singularities do not present problems (pho- 
ton energies have to be sampled as l/k°), collinear peaks, appearing in 117=1 I(h)> 
have to be carefully treated when n > 2. 

Firstly, ]T/=i 1(h) is flattened with the function Yl™ = il(h), where I(k) = 
Yh=i has the same leading singularities of the original function. We then 
write 



i=l i=l 1 

^fl^f] (— + — + — + — ) (A 5) 

1(h) jLi \Pi ■ kj p 2 -kj p 3 -kj pi-kjj 

We have to sample the second product in the right hand side of the above 
relation, which can be expanded as a sum of terms. The sum is done via Monte 
Carlo by choosing randomly for each of the photons a fermion to generate its 
angle according to 1/p- k. If rii is the number of photons "attached" to the 
fermion i, a single term of the sum takes the form 

ni i n% i 713 i n 4 -1 

l — f l — f l — f l — f / \ 



The sampling of the angular structure of the previous formula is not trivial. 
The initial-state sampling (i.e. l/pi,2 • k) is easy, being the momenta p\^ fixed, 
while to sample final-state singularities the following problem arises: consider- 
ing for example the case n = 2 with n 3 = 1 and n.4 = 1 and that we choose the 
angles of p 3 as independent variables for the phase space integral, Eq. (IA.6I) 
says we have to generate the k\ angles along the p% direction according to 
oc l/p 3 ■ ki and the k 2 angles along p 4 . The point is that the direction of p 4 
can not be defined until all the other momenta are generated. We by-pass the 
problem by using a multi-channel approach. In the case we are considering, 
we write 



1 1 1 1 pyk 

P3 -kip 4 - k 2 p 3 -fa pf k 2 -A- + 



p 3 -fci p 4 -k 2 

1 -JL_ 

1 Pi ' h2 1 (A.7) 



P3-k 1 p< l -k 2 — 1 - + 



P3-kl P4-&2 



The above sum is done again by Monte Carlo, integrating the phase space 
by using as independent variables the angles of the final-state electron (p 3 ) in 
the first term of the sum and the angles of the final-state positron (p 4 ) in the 
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second one. In the first term we generate k\ along ^3 and k% along a direction 
which is not exactly that of p±, being the l/p^ ■ k^ peak flattened in any case 
by the denominator l/p 3 • k\ + l/p^ ■ k^. In the second term the role of p% and 
P4 are exchanged. 

The generalization of Eq. flA.7|) to the case n 3 and n 4 7^ 1 reads 



1 2t 1 n 



ri3 
i=l 



tAps-k^-k, BiPs-hfApi-kj nrii^+n-ix^ 



"3 



i ft i n; 



i=l P3 • K% j=1 Pi ■ Kj H i=1 ps . fci -t- H j=1 p4 . fcj 

which suggests a similar generation of the independent variables. 
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